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^ ! Abstract 
O 

■ A microscopic theory of the transport properties of quantum point contacts 

o 

NO 

On 



giving a unified description of the normal conductor-superconductor (N-S) 



and superconductor-superconductor (S-S) cases is presented. It is based on a 

E 

i ■ model Hamiltonian describing charge transfer processes in the contact region 

and makes use of non-equilibrium Green function techniques for the calcula- 



tion of the relevant quantities. It is explicitly shown that when calculations 



^ . are performed up to infinite order in the coupling between the electrodes, 

the theory contains all known results predicted by the more usual scattering 
approach for N-S and S-S contacts. For the latter we introduce a specific 
formulation for dealing with the non-stationary transport properties. An effi- 
cient algorithm is developed for obtaining the dc and ac current components, 
which allows a detailed analysis of the different current- voltage characteristics 
for all range of parameters. We finally address the less understood small bias 
limit, for which some analytical results can be obtained within the present 
formalism. It is shown that four different physical regimes can be reached in 
this limit depending on the values of the inelastic scattering rate and the con- 
tact transmission. The behavior of the system in these regimes is discussed 
together with the conditions for their experimental observability. 



PACS numbers: 74.50.+r, 85.25.Cp, 73.20.Dx 



I. INTRODUCTION 

Electronic transport through N-S and S-S junctions and weak links has been the object 
of interest for many years Much of the theoretical understanding of these systems has 
been obtained by analyzing simple models where a one-dimensional character is assumed 
These kind of models have been very useful for clarifying the microscopic origin of 
basic phenomena like the excess current in N-S and S-S contacts and the subgap structure in 
S-S and S-Insulator-S (S-I-S) junctions. In particular, the one-dimensional scattering model 
introduced by Blonder et al. (hereafter referred to as BTK) has provided a simple way 
of analyzing the transport properties in terms of normal electron transmission and Andreev 
reflection probabilities. 

With the recent advances in the fabrication of nanoscale devices a closer com- 

parison between the predictions of simple quasi one-dimensional theories and experimental 
results is now at hand. Refs. |4]|| provide two different examples of the progress made in the 
experimental achievement of a superconducting quantum point-contact (SQPC). This pos- 
sibility has provoked a renewed interest in more detailed and quantitative analysis of models 
involving a few conducting channels. It turns out that, in spite of its apparent simplicity, 
the case of a single channel SQPC still contains non-trivial physical behavior at certain 
regimes. Single mode contacts have been analyzed recently by different authors |6|-|8[ who 
have obtained new quantitative results for the case of small applied voltage where multiple 
Andreev reflections (MAR) play a crucial role. 

Some of these recent works |7||| illustrate the existence of a variety of regimes controlled 
by the contact transmission and quasiparticle damping in the small bias limit. A full un- 
derstanding of the different regimes, together with the transitions between them remains to 
be explored. This will be thoroughly analyzed in this paper. 
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Traditionally, quantum transport in microelectronic devices has been mainly addressed 
by two different approaches: one is based on the scattering picture first introduced by 
Landauer || in which the transport properties are expressed in terms of transmission and 
reflection scattering amplitudes. This approach is generally used in a phenomenological way 
by replacing the device with a simple scattering structure. The natural extension of this 
picture to the superconducting case was provided by the already mentioned BTK model. 

A different point of view arises when the problem is analyzed starting from a microscopic 
Hamiltonian. We shall very generally call "Hamiltonian approach" the theories which take 
this starting point. The origin of this approach can be traced back to the work by Bardeen 



who introduced the tunnel Hamiltonian approximation for describing a tunnel junction [JIO 
Most of the calculations based on this tunnel Hamiltonian were restricted to the lowest 
order transport processes like in the calculation of the Josephson current in a S-I-S junction 



IIJ . Multi-particle tunneling processes (MPT) were first discussed by Schrieffer and Wilkins 
|T2|j as a possible explanation for the observed subgap structure in superconducting tunnel 
junctions. The contributions of these higher order processes were found to be divergent, 
which has led to the quite extended belief that the Hamiltonian approach is pathological 
except for describing the lowest order tunneling processes. 

One of the aims of the present work is to show that starting from a very simple model 
Hamiltonian describing a single channel contact it is possible to obtain results for N-S and 
S-S contacts in agreement with those provided by scattering theory. As will be shown below 
inclusion of higher order processes up to infinite order eliminates the pathologies associated 
with finite order perturbation theory. 

On the other hand, the Hamiltonian approach in combination with non-equilibrium 
Green function techniques presents its own advantages. In recent publications we have shown 
how this approach can be generalized for dealing with situations where self-consistency of 
the superconducting order parameter is needed [13]. Moreover, the formulation in terms of 
Green functions is specially well suited for dealing with correlation effects when strong e-e 
interactions are present ||14|| . 



The paper is organized as follows: In section II the general model for a superconducting 
weak link in a local representation is introduced. In particular, we discuss how this model 
can be simplified in order to describe a single mode contact. We then introduce the non- 
equilibrium Green function formalism by means of which the total current through the 
contact can be expressed. In section III we analyze the more simple N-N and N-S contacts. 
The study of the N-N case allows us to define the contact transmission coefficient in terms 
of the microscopic parameters of our model. The complete correspondence with the results 
of the BTK scattering theory is then established by analyzing the N-S contact. For this 
last case we derive an analytical expression for the total current as a function of the contact 
transmission. Section IV is devoted to the S-S case. We first show how the problem of the 
calculation of the ac current components can be reduced to the evaluation of an algebraic 
set of linear equations. This allows us to study in detail the general features of the dc 
and ac / — V characteristics for the whole range of voltages and transmissions. We finally 
concentrate on the limit 7 < A where the more interesting and less understood physics 
takes place. Section V is devoted to some concluding remarks. 



II. MODEL AND METHOD 

In this paper we shall consider "point contact-like" geometries consisting very generally 
of two wide (3D or 2D) electrodes connected by a narrow constriction. The constriction 
length, Lq, is assumed to be much smaller than the superconducting coherence length and 
its width, Wc, comparable to the Fermi wavelength , Xp. For a sufficiently short constriction 
the detailed form of the self-consistent order parameter and electrostatic potential in the 



region between the electrodes becomes irrelevant ||T3| allowing us to represent them by 
simple step functions. On the other hand, the condition Wc ~ Xf implies that there is 
only a reduced number of quantum transverse channels through the constriction. One can 
think of different experimental realizations of such a physical situation, like for instance the 
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recently developed atomic size break-j unctions f|] and the split-gate S-2DEG-S SQPC of 
Takayanagi et al. |J. These two situations are schematically represented in Fig. 1. 

The analysis of the transport and electronic properties of such a system would involve, 
for the most general case where both electrodes are superconductors, the solution of the 
time-dependent Bogoliubov-de Gennes (BdeG) equations fL5| . In previous works we have 
developed an approach to these kind of problems in which the BdeG equations are formulated 
in a site representation |fT3| . This representation can be viewed either as a tight-binding 
description of the electronic states or as a discretization of the BdeG equations. The first case 
would be more suitable for describing systems like an atomic size contact (Fig. 1 a), while 
the second could be used to represent constrictions involving a 2DEG in semiconducting 
heterostructures (Fig. lb). In the absence of applied fields the mean- field Hamiltonian 
giving rise to the BdeG equations would take the form 

h = XX e * - mKW + H tijcUj + H( a * c !i c !t + A iC it c a ), (i) 

where i,j run over the sites used to represent the system. By appropriately choosing the 
different parameters e,, and Aj one can model different junction geometries [13]. In order 



to analyze the case of biased contacts it is convenient to simplify this model Hamiltonian 
taking the step-like behavior of the order parameter and electrostatic potentials in a point 
contact geometry into account. To do this we take the complex order parameter and the 
chemical potentials as constants on the left and right electrodes (denoted by (Al,/!^) and 
(A R ,fi R ) respectively). The case of a single quantum channel connecting both electrodes 
can then be described by the following Hamiltonian 

H = H L + H R + J2(tcLc R * + t*c Ra c La ) - /2 L N L - fj, R N R , (2) 

cr 

where Hl and H R correspond to the uncoupled electrodes while the hopping term describes 
the electron transfer processes between the outermost sites on both electrodes. We would 
like to stress that, although this is a very simple model Hamiltonian (formally equivalent to 
a tunnel Hamiltonian), it contains the relevant physics of a quantum point contact which 



depends essentially on the contact normal transmission coefficient. In our model this trans- 
mission can vary between ~ (tunnel limit) and ~ 1 (ballistic regime) as a function of the 
coupling parameter t. This will be discussed in section III. 

For the case of a symmetric superconducting contact in the presence of an applied bias 
voltage eV = — //r it is convenient to perform a gauge transformation [l6j by means of 
which Hamiltonian (2) adopts the following time-dependent form 

H{r) = H L + H R + Y: {te^c\ a c Ra + t*e~^% a c Lf7 ) , (3) 

where the superconducting phase difference (j>{r) = <fio + 2eVr/h appears as a phase factor 
in the hopping elements. 



A. Current in terms of the non-equilibrium Green functions 



Within the previously introduced model, the average total current through the contact 
is given by 



Hr) = ~ £(* < cUt)crAt) >-t*< cL(r)c L(T (r) >) 



(4) 



where, depending on the gauge choice, t can include a time dependent phase like in Eq. (3). 
The averaged quantities in Eq. (4) can be expressed in terms of non-equilibrium Keldysh 



Green functions G + | 17[| . For the superconducting state it is convenient to introduce the 



2x2 Nambu representation in which G]~ 3 - adopts the form |T3[| 



< cjt( r, ht( r ) > < 9 i (r')Q T (r) > 

^ < 4t( t/ ) c !i(^) > < ^(rOd^r) > ) 

In terms of the G + ~ , the current is given by 

2e 
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(5) 



(6) 



where t is the Nambu representation of the hopping elements 



t 



) 



(7) 



t 

-t* 

The problem then consists in the determination of these non-equilibrium Green functions 
G + ~ . This can be formally achieved by treating the coupling t as a perturbation. The 
unperturbed Green functions correspond to the uncoupled electrodes in equilibrium. For 
a symmetric contact, and neglecting finite band-width effects, the uncoupled retarded and 
advanced Green functions can be expressed as 



9llM = 9rr{u) 



\ 



V 



-oj ± ir\ A 
A —ijj ±ir] j 



(8) 



W^JA 2 -(u± irj) 2 

where W is an energy scale related to the normal density of states at the Fermi level by 
p(€f) ~ 1/ (ttW) and 77 is a small energy relaxation rate that takes into account the damping 
of the quasi-particle states due to inelastic processes inside the electrodes [ITH] . On the other 
hand, the unperturbed g + ~(uj) satisfy the relation 

g + ~{oj) = 2nip(u)n F (u), (9) 

where p(uj) = (l/iT)Im[g a (uj)} and np{uj) is the Fermi function. 

For the coupled system, the functions G^J can be obtained from the retarded and ad- 
vanced functions by means of the integral equation 



G+-(t, r') = J dndT 2 [I5(r - r x ) + G r (r, n^in) 
'i5(r 2 -r') + ± a (r 2 )G a (T 2 ,r') 



9^ [n - t 2 ) x 



(10) 



where G r,a can be derived from their corresponding Dyson equations 

G r > a (r, r') = g r ' a (r - r') + J dng r ' a {r - T 1 )tr> a (T 1 )Gr> a (n, r') 



(11) 



In the above equations the self-energies take the simple form S£ L = = and 



In sections III and IV details of the solution of these integral equations are given. 
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III. N-N AND N-S CONTACTS 



In this section we briefly review the N-N and N-S cases. The analysis of the N-N case 
permits to define the normal transmission coefficient of the contact in terms of the micro- 
scopic parameters of the model, which is necessary for making contact with the scattering 
approach. The equivalence of both approaches is illustrated by comparing our results for 
the N-S contact with those of the BTK model @. 

As we shall deal with a stationary situation it is convenient to adopt the time-independent 
formulation based on Hamiltonian (2). In this case the Green functions depend on the 
difference of its temporal arguments and integral equations (10) and (11) become simple 
algebraic equations when Fourier transformed. Then we have 



G r ' a {u) = g r ' a {uj) + g r ' a (cu)± r ' a (uj)G r ' a (iu) 



(12) 



G< + -)'(- + >(w) = [l + G r (u)£ r {uj) 
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(+-),(-+)/ 



u) I + t a {u:)G a {uo) 



(13) 



where g = — 2-nipiiS) [1 — np(uj)]. On the other hand, it is convenient to rewrite 

expression (6) for the current as [ |l3| . |Kf] 



2e 9 ' N 



9LL,\\G RRAli^) ~ 9 'll, 11 C WnO^O 



(14) 



where we have assumed that the left electrode is in the normal state, while the right one 
can be either normal or superconducting. 

Let us first analyze the case where both electrodes are in the normal state. Details on 



the calculations of the Green functions for this case have been given elsewhere [ 2~0fl . As in 
the scattering approach, the current can be written as 

2e roc 



ze f°° 

I — — / T(u, V) \jif\u — eV) — n F (uj)\ du, 

tl J —oo 



(15) 



where T(u, V) is an energy- dependent transmission coefficient which is given by 

4?r 2 t 2 p LL (to - eV)p RR (u) 



T(u,V) 



\l-Vg LL (uj-eV)g RR {uj)f 



(16) 



This expression can be further simplified assuming that the normal state spectral densi- 
ties are constant over an energy scale W 3> V, A. This is consistent with the assumptions 
leading to Eq. (8). Within this energy scale the coefficient T(u>,V) becomes an energy- 
independent quantity a 

At 2 IW 2 

r <^(TT7W¥^ (17) 

The normal linear conductance of this single mode contact is then given by the Landauer 
formula G^n = (2e 2 /h)a. Notice that a can vary between zero and one as a function of t. 
The a — > limit is reached both for t/W <C 1 and for t/W 3> 1, while the ballistic limit, 
i.e. a ~ 1, is reached when t/W ~ 1. 

It is worth to further clarify the role of the hybridization parameter t within our point 
contact model. Although in the context of the tunnel Hamiltonian approach it has been 
customary to identify ~ t 2 with the transmission probability, this strictly holds for the 
tunnel regime (lowest order perturbation theory in t). However, the actual expression for 
the contact transparency (Eq. (17)) including all order processes is a non-linear function of 
t 2 . While the ballistic condition is achieved for t 2 /W 2 ~ 1, the way in which a approaches 
unity is completely different from ~ t 2 /W 2 , as can be seen from Eq. (17). 

In the N-S case (A L = 0, A R = A) starting from Eq. (12) and (13) after some sim- 
ple algebra we obtain the following expression for the current as the sum of four different 
contributions I — I± + I 2 + I 3 + Ia, where 
8e r°° 

h = -r^H 2 / du\l + tG r RL ll {uj)\ 2 p LL>11 (u - eV)p RRM (u) [n F (uj - eV) - n F (uj)} 
ri J —oo 

16c f°° r r i i 

h = y_ oo d " Re {tG a LR , 21 (u;) [l + tG RL , n (u)} } x 

Pll,\\{u - eV)p RRtl2 (uj) [n F (uj - eV) - n F (u)} 

gg rco 

h = -r^H 4 / du\G RLil2 (u)\ 2 p LL!l i(uj - eV)p RRt22 (uj) [n F (uj - eV) - n F (u)\ 
ri J -co 

gg roc 

I A = -t-ttH / du;\G RR)12 (u)\ 2 p LLtll (uj - eV)p LLt22 (u + eV) [n F (uj - eV) - n F (u + eV)] . 
ri J -co 

(18) 

Written in this form, each contribution has a clear interpretation in terms of elementary 



processes that can be identified by inspection of the intervening spectral densities. Although 
there is a formal resemblance of the above expressions with those of tunnel theory JZT] Eq. 
(18) contains all possible processes up to infinite order in t. Thus, I\ corresponds to normal 
electron transfer between the electrodes, while I2 corresponds also to a net transfer of a 
single electron with creation or annihilation of pairs as an intermediate state. On the 
other hand, ^3 arises from processes where an electron in the normal electrode is converted 
into a hole in the superconducting side, i.e. processes with "branch crossing" in the BTK 
language. Finally, I a arises from Andreev reflection processes in which an electron (with an 
associated spectral weight pll,ix{w — eV)) is transmited from the left to the right electrode 
with a hole reflecting backwards into the normal electrode (with an associated spectral 
weight pLL^iyi + e ^0) while a Cooper pair is created in the superconducting side with a 
probability proportional to {Grr^uj) | 2 . 

As expected, the only non-zero contribution for eV < A is Ia, while all processes con- 
tribute for eV > A. With the same simplifying assumptions leading to Eq. (17) the 
differential conductance at zero temperature adopts the simple form 

4e 2 a 2 

G NS (V) = ^- — ; eV > A. (19) 

h a + (2-a) v /l-(A)2 

This expression can be shown to be equivalent to the one obtained from the BTK model 
with the correspondence Z = [1 — (t/W) 2 ]/(2t/W), where Z is the dimensionless phe- 
nomenological parameter controlling the barrier height in the BTK model. 

It is worth noticing that the differential conductance rises from (4e 2 /h)a 2 / (2 — a) 2 at V = 
to the value 4e 2 /h at eV = A, this last value being independent of the contact transmission. 
This result can never be obtained within any finite order perturbative approximation in t, 
but requires an infinite order calculation [22 ] . This sort of non-perturbative features are also 



very important in the S-S case [E3J, as will be discussed in the next section. 



Finally, it is possible to obtain from Eq. (18) an analytical expression of the total current 
at zero temperature. For eV > A this can be written as / = I\ + 12 with 



in 



h 



h(x) 



eA 



a" 



h {2-a)^/T^a 



In 



1 + 






2-Q 


1 - 






2-a 



4e 

77 



A 



4(1 -a 



+ 



a 



a(2 — arx 



x 



a 



+ (2-a)vT 



a 2 (2-a) 
+ 8(1 - a) 3 / 2 



l+v^ 



+ 



1 



4(1 -a) a + {2-a)Vl~x 2 



(20) 



where ii is the gap contribution to the total current, I2 is the contribution coming from 
energies outside the gap and x = A/eV. This allows to analyze with detail the "excess" 
current, defined as I exc = limv-_ >0O (IjV'S — Inn), as a function of the contact transparency. 

~2JT^a~ 



We find I exc = I exci 



eA 



a 



h (2-a)y/T^ r a 



In 



1 



+ 



1 



a 



2-a 



2J\^x 



2-a 



In 



1 - VT^a 
1 + VI - a 



I- a 2(1 -a) 3 / 2 

where I exCl and I exC2 are respectively the contributions coming from energies inside and 
outside the gap. As can be easily checked, I exci > while I exC2 < the total excess current 



being always positive. For a = 1 the well known result [24, 28 1 for the ballistic contact 



Iexc = (8/3)eA//i is recovered. 



IV. THE S-S CONTACT 



A. An efficient algorithm for evaluating the ac current 



As commented in section II, for the case of a voltage-biased S-S contact it is convenient 
to start from Hamiltonian (3) in which the applied bias is taken into account through a 
time-dependent phase factor in the hopping element, which in a Nambu representation has 
the form 



t 



te iHr)/2 

-t e ^ (r)/2 



(21) 



1 1 



where = <fto + 2eVr/h is the time-dependent superconducting phase difference. 

This explicit time dependence indicates that all dynamic quantities can be expanded as 
Fourier series in all possible harmonics of the fundamental frequency o;o = 2eV/h J7|. For 
instance, the total current can be written as 

I(T) = Y,I m e muJ0T . (22) 

m 

We shall now show how these Fourier coefficients, I m , can be efficiently evaluated within 
the non-equilibrium Green function formalism. Let us first notice that the non-equilibrium 
Green functions appearing in Eqs. (10) and (11) do not depend only on the difference of 
their temporal arguments, and have therefore a generalized Fourier expansion of the form 



25. m 



G(t, t') = — V / duje- luJT e i{uJ+nuJo/2)T 'G{uj, u + ntv /2). (23) 
2n n J 

Hereafter we shall use the notation G nm (uj) = G(u>+nu /2, u+muj /2). Different Fourier 
components G nm are related by G nm (uj) = G n _ m>0 (^ + mu> /2). For the particular gauge 
choice adopted here, it is useful to express all quantities in terms of a renormalized hopping 
which satisfies its own Dyson equation 

f a > r (r, rO = t(r)5(T - r') + J dndr^g^ij - n)& (n - r 2 )f^{r 2 , r'). (24) 

This quantity can be viewed as the total hopping amplitude arising from summing up 
all processes in which one electron is transferred. Clearly, it is formally equivalent to use a 
renormalized hopping instead of renormalized propagators as they are linked by relations like 
GLL( T , T ')i( T ') = I ^ t \9l[j — t \)Tlr{ti,t')- The current components can now be expressed 
in terms of the renormalized hopping elements T^(u) as 

T - — fr/wV [f r o+-f ,r t rt a _ P,T fr ~+-fri 
h I ' ' I OlffiTl ± nmiimm yOO^Oni/nn ± nm 

II J n 

where we have eliminated the site indexes L and R in the uncoupled Green functions due 
to the left-right symmetry of the contact. 

1 9 



The problem is then reduced to that of the evaluation of the components f nm . From Eq. 
(24) it can be shown that the components T nm (both retarded and advanced parts) satisfy 
a set of linear equations of the form 

Tnm tnm ~\~ ^nXnm ~l~ ^4,(1-2^1-2,171 ~l~ ^ / n,n+2^~n+2,m- (26) 

These equations are mathematically equivalent to those describing the motion of electrons 
in a tight-binding linear chain with "site energies" , e n , and "nearest-neighbor couplings", 
V n ,n-2 and V n ,n+2- The detailed expression of e n and V n) m in terms of the unperturbed Green 
functions are given in Appendix A. This analogy allows us to obtain the Fourier coefficients 
T nm using standard recursion techniques. One can show (see Appendix A) that the following 
recursive relation holds 

T n +2, m (u) = z + [u + (n - l)u ] f nm (u) , n > 1 

v^7) 

T n _ 2: m{oj) = z~ [u + (n + l)u ] T nm {uj) , n < -1 
where the transfer matrix satisfy the equation 



-i 



(to) =[I- e ±3 - V± 3 ,±5^ (oj ± uj ) \ ■ (28) 

Clearly, as the transfer matrix connects consecutive harmonics of T, it can be viewed 
as a generating function which introduces the effect of a unitary Andreev reflection process. 
The problem has been reduced to the calculation of only two matrix coefficients like, for 
instance, T lj0 and T_i )0 as a starting point for the generating Eqs. (27) (see Appendix A for 
details). 

In summary, the basic mathematical difficulty lies in the evaluation of the transfer matrix 
functions z^ from Eq. (28). Although Eq. (28) looks simple, it is nevertheless hard to solve 
analytically for arbitrary values of V . The analytical results presented in this paper are 
limited to the eV/A — > and eV/A — > oo cases where some simplifying relations hold. For 
intermediate voltages, an accurate numerical solution of Eq. (28) can be obtained. 
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B. Analysis of ac and dc I — V characteristics 



In this subsection we analyze the general features of the I — V characteristics of a S-S 
contact obtained using our formalism. Let us start by briefly discussing the dc current, 
Io{V), for different values of the transmission, as shown in Fig. 2. Although the overall 
qualitative features of these curves have been known since the works of Octavio et al. p7| . 



Zaitsev [24] and Arnold p5| , more quantitative and detailed analysis are being reported in 
recent publications . The results of Fig. 2 are in agreement with those reported recently 
by Averin and Bardas || which were obtained using the scattering approach. 

Two relevant features of these curves are the subharmonic gap structure for eV < 2A 
and the excess current for eV ^> A. As can be observed, the subgap structure becomes 
progressively more pronounced with decreasing transmission. Eventually, when aCl the 
current steps at positions eV ~ 2 A/to can be clearly resolved. In this limit one can isolate 
the elementary processes which give rise to these steps. It can be shown that in the tunnel 
limit the n-th step can be calculated as 



S P r neV-A 



Sir = lim — ir 2 nt 2n / du 

eV^2A+/n h J A 



n-1 

pn{oS)p^{u — neV). (29) 



H\g 12 (u-ieV)f 



.i=i 

By comparison of Eq. (29) with the expression of Ia in Eq. (18) it becomes clear that 
the steps inside the gap are due to the opening of a new Andreev reflection channel whenever 
eV = 2 A/to. Calculation of the integrals in Eq. (29) leads to 

(n) eA«" / 2to \ (n n \ 2 

in agreement with the recent prediction of Bratus et al. 0. On the other hand, when a — *■ 1 
the subgap structure is completely washed out and there appears an excess current even in 
the small bias limit. 

The other relevant feature of the dc current, namely the large voltage excess current, 
can be analytically evaluated within our model for any transmission value. The main sim- 
plification in this limit comes from the fact that only the lowest order Andreev reflection 



i/i 



process gives a significant contribution to the excess current ||27|| . This implies that one can 
truncate the system of equations (26) for harmonic indexes n > 1, the resulting simplified 
system can then be solved explicitly for T 10 (co>) and T„ 10 (c<;) (for details see Appendix B). 
As shown in this Appendix the simple result l££ c = 2/^f is obtained for any value of the 
transmission. Although this is the expected physically sound result, to our knowledge, it 
has not been shown explicitly before except for the ballistic case p4p8f . Moreover, some 



authors have reported the existence of a negative excess current for low transmissions f29 
which seems to be in contradiction with the above result. However, one should notice that 
the excess current as defined above is an asymptotic quantity (only valid in the eV/ A — > oo 
limit). When corrections of order A/eV are taken into account one actually can have a 
defect instead of an excess current for sufficiently low transmission. 

The algorithm described in subsection A allows an efficient evaluation of the higher order 
ac components of the current. For the following analysis we decompose the ac current, Eq. 
(22), into its dissipative and non-dissipative contributions given respectively by 

I D = h + Y.C^{muj Q T) (31) 

m 

and 

J 5 = £Z£sin(ma;oT), (32) 

m 

where 1% = 2Re(I m ) and 7^ = -2Im(I m ). 

The results obtained for the first three I® and 1^ components are depicted in Fig. 3 and 
Fig. 4. As can be observed, these components become exponentially small for bias voltages 
larger than A/n. On the contrary, when eV < A/n the decay of the ac components with 
increasing n becomes slower. The analysis of the higher order components reveals a decay 
for eV < A/n close to a inverse power law. As a consequence of this slow decay, one is 
forced to take an increasing number of ac components into account in order to adequately 
describe the behavior at small bias. This will be the subject of the next subsection. 



C. Small bias regime 



In this subsection we concentrate on the eV/A — > case, which turns out to exhibit a 
remarkable variety of different regimes according to the values of the parameters qA and the 
inelastic scattering rate rj. As it is well known, the main difficulty for obtaining quantitative 
results in this limit lies in the fact that the number of MAR contributing to the current grows 
with decreasing V as ~ A/eV |p5| . Furthermore, the amplitudes of these multiple processes 
do not decay when (V, rj) — > leading to the appearance of divergencies in the perturbative 
expansion in the coupling t 0. Again, a complete summation of the perturbative series is 
needed in order to regularize these divergencies. An additional difficulty arises, as will be 
discussed below, from the fact that the limits V — > and rj —>■ do not actually commute. 

When decreasing the bias voltage two main situations can be reached depending on 
the strength of the inelastic scattering rate at the superconducting electrodes: the case of 
eV I A — > and finite rj has been discussed by the present authors in recent publications [|7| ; 
on the other hand, the case of small V/ A and negligible rj has been recently addressed to 
by Averin and Bardas |§. In what follows we summarize the main results for both regimes 
and analyze the conditions for their actual observability in a real SQPC. 

Within our formalism, analytical results in the small bias limit become feasible since the 
transfer matrix z (oS) tends to a scalar quantity having the form of a simple phase factor 
inside the gap region. As shown in Appendix C for (77, V) — > we find 

z ± (u) = z(w) = e iip(u) , Ay/T^ < M < A, (33) 

where 



ip(u) = arcsin ( — — a/A 2 — uj 2 \ uj 2 — (1 — a) A' 
\aA 2 v 

As the multiple Andreev processes are generated by successive applications of z(u), Eq. 
(33) indicates that these processes do not decay in this limit inside the gap region. This 
infinite series of MAR gives rise to the well known bound states spectrum of a current- 



carrying SQPC at zero bias voltage |p0| . As shown in Appendix C, the positions of these 

1 fi 



bound states are determined by the condition <p(uj) = <fi. The presence of a small but finite 
1] or V (whichever is larger) introduces an effective damping into the otherwise infinite series 
of MAR. 

When this effective damping is due to a finite rj (eV <C 77) a linear regime can be defined 
where the total current is given by Is{4>) + Cf(<j>)V, G(4>) being a phase-dependent linear 
conductance |]24lj7|1. Within this linear regime, one can identify two different sub-regimes 
according to the ratio rj/aA. The case rj/aA <C 1 corresponds to a situation where MAR 
are very weakly damped and give the dominant contribution to the current. The physical 
picture one can have in mind is that of electrons and holes Andreev-reflecting between the 
electrodes for a very long time before being inelastically scattered |31| . In order to illustrate 



the dominant contribution of the processes inside the gap for the weakly damped case, we 
represent in Fig. 5 the current density corresponding to the Iq component for three values 
of the transmission. Three important features of this linear regime are displayed in this 
figure: firstly, the current density inside the gap increases as ~ l/i] therefore giving the 
dominant contribution in the weakly damped case; secondly, there is a region inside the 



gap of width 2Aa/1 — a in which the current density vanishes. This is the forbidden energy 
region for bound states at a given transmission. Finally, in Fig. 5 (c) one can observe 
that the contribution of the continuum outside the gap becomes important as a < rj/A. In 
the limit r]/aA > 1 a second sub-regime is reached where the contributions of MAR are 
heavily damped and the current is dominated by single quasi-particle tunneling processes. 
The transition between these two sub-regimes has been analyzed with detail in Refs. 0. 
In order to identify the actual sub-regime for a real SQPC an estimation of the order of 



magnitude of r] is needed. In Ref. [|18[] rj is estimated from the electron-phonon interaction 
to be a small fraction of the gap for traditional superconductors. Thus, our theory predicts 
that a real SQPC would generally fall into the weakly damped case except for extremely 
low transmissions. 

For this sub- regime the supercurrent Is and the linear conductance G(4>) can be obtained 
analytically as discussed in Appendix C. In particular, for G((p) one obtains 

1 7 



CM 3 ' : ' 7 



h 16r] 



Aa sm<f) ,(3ojs* 
=sech( — - 
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P, (34) 



-V /l - a sin 2 (0/2) ' 

where u>s is the position of the bound states inside the gap and (3 = l/ksT. This expression 
for large transmission and small temperatures gives a phase-dependence which is in quali- 
tative agreement with the few available experimental results, performed in non-mesoscopic 
contacts The unusual phase- dependence of Eq. (34) which deviates strongly from 

the cos(0) form predicted by the standard tunnel theory may explain the old controversy 
between tunnel theory and experiments known as the cos(0) problem [0,0 • 

On the other hand, when the truncation of the infinite series of MAR is caused by a 
finite V (with negligible 77), analytical results have only been obtained in the quasi-ballistic 
limit, i.e. a — > 1. A closer inspection of the I — V curves in the small bias region and 
for a ~ 1 reveals that the supercurrent components decay exponentially from its value at 

V = with a collapsing width ~ (1 — a) A. This is illustrated in Fig. 6 (a) where a blow up 
of the behavior of If for small bias voltages and quasi-ballistic transmissions is shown. In 
the limit a — > 1 the supercurrent becomes a delta function at V = 0. On the contrary, the 
dissipative components in this same limit tend to a finite value outside the region of width 
~ (1 — a) A. This behavior is shown in Fig. 6 (b) where if is plotted in the same magnified 
scale as If. The summation of these dissipative components for a = 1 and very small V 
yields (see Appendix C) 

eA 

I D (<f ) ) = —\sm( ( j ) /2)\signV, (35) 

in agreement with the result recently derived by Averin and Bardas ||. The existence of 
a region of decreasing width V ~ (1 — a) A in which this crossover from supercurrent to 
dissipative current takes place can be associated with the collapse of the forbidden region 
for MAR inside the superconducting gap taking place when a — > 1. In this way, when 

V is small compared to the width of the forbidden region the excitation of quasi-particles 
from states at uj < — A^l — a into states at uj > A^l — a is negligible and there is no 
appreciable dissipative current; whereas the opposite situation holds for V > y/l — «A. 

is 



Averin and Bardas || have described this crossover as a Landau-Zener transition in which 
the non-dissipative and the dissipative components scale with a and V as (1 — p) and p 
respectively, where p = exp [— 7r(l — a)A/eV]. The numerical results for sufficiently small 
eV/A and (1 — a) are well fitted by these scaling laws. However, a careful analysis reveals 
that their range of validity around V = and a = 1 decreases strongly when increasing the 
component number. 

In summary, in this small bias limit one can identify four different sub-regimes depending 
on the relative values of parameters r], aA and eV. The prediction of the actual behavior of 
a real SQPC in this limit would therefore require a careful estimation of all these parameters. 
In this respect, one should keep in mind that while eV and a can be varied experimentally 
in a rather controlled way, the inelastic scattering rate 77 is an intrinsic property of the 
superconducting electrodes much more difficult to control. The unavoidable presence of 
some degree of inelastic scattering can prevent the actual observability of the crossover 
from non-dissipative to dissipative behavior described after Eq. (35). The requirement of 
eV ~ (1 — a)A together with that of a ~ 1 can actually imply eV < rj which would rather 
correspond to the linear regime. Finally, when considering the experimental test of all these 
theoretical predictions the relevance of noise in a real SQPC should be taken into account. 
Recent theoretical predictions f35|,|36| suggest that the magnitude of thermal noise in a single 
mode superconducting device can be extremely large near ballistic conditions. 



V. CONCLUDING REMARKS 

In the present work we have presented a Hamiltonian approach for describing the trans- 
port properties of single mode N-S and S-S contacts. 

It has been explicitly demonstrated that this approach is, with some simplifying as- 
sumptions, equivalent to the phenomenological scattering approach. We believe that the 
present work can help clarifying the somewhat recurrent discussion about the unsuitability 
of a Hamiltonian approach for obtaining the transport properties of a N-S or S-S contact: 



1 



when performing the calculations up to infinite order in the coupling t all the unphysical 
divergences are eliminated and the results become equivalent to those of scattering theory. 

On the other hand, the present approach has been applied to discuss in detail the dc and 
ac I—V characteristics of a SQPC. In particular, we have concentrated on the less understood 
small bias voltage limit where one can identify four different sub-regimes depending on the 
values of the contact transmission and inelastic scattering rate. Finally, we have discussed 
the conditions for the experimental observability of the theoretical predictions in this limit. 

Although in the present work we have restricted the discussion to the simplest single- 
mode case with an energy- independent transmission coefficient, the general model intro- 
duced in section II can describe more complex situations which may be relevant for a closer 
comparison with recent experiments || . Work along this line is under progress. 
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APPENDIX A: 



In this Appendix we describe the algorithm for evaluating the ac current components in 
a S-S biased contact. As it has been pointed out in section IV, we can express the current in 
terms of the retarded and advanced Fourier components of the renormalized hopping T nm {uj) 
satisfying Eq. (26) 



T n m trim ^ri^nm *n,n— 2,m K!,n+2^n+2,m) 



(Al) 



where the matrix coefficients e n and V nm can be expressed in terms of the uncoupled Green 
functions as 



9n+l9n 9n+l fn 
9n-lfn 9n-\9n 



(A2) 
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Vn,n+2 — —t 2 fn+l 



Jn+2 9n+2 



V 



o o 



(A3) 
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Ki,n-2 — ~t 2 fn-l 



(A4) 





^ 9n-2 fn-2 j 

where f{uj) = #12 (k>) = 5 , 2i(<^) and g{uj) = gu(u>) = #22 (u>) due to electron- hole symmetry. 
In the above equations the short-hand notation g n = g{uj + noo /2) is used. Moreover, the 
site indexes in the Green functions have been omitted since we are considering a symmetric 
contact. 

As commented in section IV, the linear Eqs. (Al) are analogous to those describing a 
tight-binding chain with nearest-neighbor hopping parameters V n ,n+2 and V n ,n-i- A solution 
can then be obtained by standard recursive techniques. It is straightforward to show that 
the following recursive relations between the coefficients T nm hold 



T n+2 , m (u) = z + [uj + (n - l)u ] T nm (uj), n > 1 
f n - 2 ,m(u) = z~[uj + (n+ l)u ] f nm (uj), n < -1 

where the transfer matrix ^(cu) satisfy the equation 



(A5) 



z ± (u) = I - e ±3 - V ±3t ± 5 z ± (uj±u} ) 



i -1 



(A6) 



One can see from Eq. (A6) that z + {uj) and z~{uj) are related by z~(u, V) = & x z~(u, —V)a x , 
where a x is the corresponding Pauli matrix. 

By virtue of the relation T nm {uj) = T n „ mi0 (^ + mcJ /2), one can write the current com- 
ponents given by Eq. (25) in terms of T n0 (uj) = T n . Using recursive relations (A5) the 
calculation can be reduced to a closed system for coefficients T\ and T_i 



/ - ei - V 13 z + (cu)\ Ti = £10 + 
/ - e_i - V- lt - 3 z- (cu)} f_! = L 10 + VLliTl 



(A7) 
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The remaining task is the calculation of the transfer matrix z + {uj). It can be shown that 
the solution of Eq. (A6) is a diagonal matrix whose elements can be expressed in terms of 
a scalar function \ + {uj) 

\ 



z+(u) = -t 2 



( f f s o s t 







V 



5+ 8+ 

'■'■#/ 



(A8) 



where A+ = + nu /2) and 5+ = (A+ — g n +2)/t 2 f 2 +2 - The function A + (cj) satisfies the 
following equation 



A+ = a + b\£ + c\t + d\$\$, 



(A9) 



where, taking into account Eq. (17), the coefficients can be written as a = g 2 + (t 2 /W 2 )g 3 , 
b = t 2 g 2 g 3 , c = — t 2 [g 2 g 3 — (t 2 /W 4 )] and d = t 2 [g 3 + (t 2 /W 2 )g 2 ]. For an arbitrary bias voltage, 
Eq. (A9) can only be solved numerically. However, as we show in Appendix B and C it is 
possible to obtain analytical solutions in the two special cases ujq — > and u — > oo. 

Once z + (uo) has been determined one can calculate the coefficients T\ and T_i from Eq. 
(A7) 



Ti = 



-f 



1 t'hhb 2 b\ t 2 hb\\ 



\-t 2 \ 2 \\ 



t 2 fo5 2 



(A10) 



f. 1 (u,V) = -a x f 1 (u,-V)a x , 



(AH) 



where \ n (uj,V) = A+(w, — V) and 5 n (u, V) = 8+{u,—V). The rest of the coefficients f n 
can be calculated from Eqs. (A5) 



- 2n+l 



Li=l 



n > 



f-n = -a x f n (uj, -V)a x , n > 0. 



(A12) 



Finally, the current components, separated into its dissipative and non-dissipative parts 
(Eqs. (31) and (32)) can be calculated from the expressions 
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4g roc 

h = - T ] dw £ Re{Tr(a z Tl(uj)^-T n (u)g a )} 



(A13) 



n=odd 



I ™ = ~JiLoo duJ ^ Re{Tr(a z [fl +m (u-mu /2)+fl_ m (u + rruv /2)]g+ T n (cc>)(?o)} 



n=odd 



(A14) 
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n=odd 
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APPENDIX B: 

In this Appendix we give details on the evaluation of the excess current for S-S contacts. 
In the limit eV/A — > oo only the dc current component Jo survives. The infinite summation 
over n in Eq. (A13) can be truncated in this case neglecting the n > 1 terms. This is 
justified as the products f n fn+i are negligible in this limit leading to a vanishing transfer 
matrix £ ± (o;). Physically, this is equivalent to neglecting multiple Andreev processes for 
eV/A > 1. Then, Eq. (A13) reduces to 

4e r°° 

du, 

71= -1,1 



h 



du £ Re{Tr(a z fl(u) 9 :-(u)f n (u;)g a )}, 



(Bl) 



with 



71 



' t 2 /i^±i ^ 



1 - t 2 \ 2 \t 1 



t 2 f 5 2 



T^(uj,V) = —a x Ti(u, —V)a x . 

On the other hand, when neglecting contributions of order A/eV Eq. (A9) simply yields 
K ~ (#71+2 + it 2 /W 3 )/(l - it 2 g n+2 /W). We then obtain from Eq. (Bl) the simple result 
ifxc = ^exc -i f° r the excess current at zero temperature and any value of the transmission 
coefficient. 
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APPENDIX C: 



In this Appendix we give the main steps in the analytical calculation of the current 
components in the limit eV/A — > 0. 



Linear regime (77 3> eV) 

The small voltage response can be straightforwardly derived from Eqs. (A13-A15) by 
expanding the Fermi functions appearing in g + ~ up to first order in eV: npiuj + nu /2) ~ 
n F (uj) — ([3 /d>)nuj secYi 2 (f3uo/2) and evaluating the rest of these expressions at eV = 0. The 
current components can be then written as 

2e 2 _ r°° , , 



Io = ^~(3V r Aosech 2 ^) £ nRe {Tr(a z fl(g a - g r )f n g a )\ 



(CI) 



n=odd>0 



2e 2 



C=~-J-PV Aosech 2 ^) £ ni?e{Tr(a 2 

n=odd>0 



7-4 4. 'ft 
^-n+m 1 - 1 n- 



(g a - g r )f n9 a )} (C2) 



gg /-oo . 

I " l = ~h J -00 du ^ /m l Tr ^ 2 



'ft 'ft 

- 1 n+m ± n—m 



(g a - g r )f n g a )} . (C3) 



8e f 00 

n=od(i>0 

Thus, the dissipative contribution, 7 D , goes to zero as Id{4>) ~ being the 

phase-dependent linear conductance, while the supercurrent part, 15(0), tends to a finite 
value at V = 0. 

In the zero voltage limit the coefficients T n adopt a simple form. The transfer matrix 
^{uj) becomes a scalar function: z + {uj) = z~{uj) = z(u)I; with = —t 2 f5 2 /X 2 , where 
A + (c<j) = A~(cj) = A (a;) satisfies the simple quadratic equation 



t 2 



^HA 2 H-(i-— )ah + j( w ) = o. 



(C4) 



Finally, the coefficients T„ adopt the form 

-t 



1 - t 2 A 2 



( t 4 f 2 5 2 t 2 f5 
t 2 fS 1 



(C5) 
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T 2n+1 (cu) = z n (cj)T 1 (cu), n > 0. 



(C6) 



Due to these simple recursive relations, the series appearing in the current components 
become geometrical series, which can be summed up without difficulty. In the weakly 
damped case, rj/aA <C 1, these summations lead to analytical expressions for the dissipative 
and non-dissipative parts of the current. By solving Eq. (C4) up to corrections of order 
rj/aA one obtains 



z{uj) = e 



Acut) r . 



«A 2 



i + cotg(<p(w))] , AVI - a < \lo\ < A, 



where 



<p(u) 



arcsm 



The summation of the geometrical series yield 
2e 2 „„ r 00 , , ,Pu 



da; sech z (^— )i?e {^(a;)} 

oo 2 



# = du> sech 2 (^)i?e | 
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where = Tr[a z T;(g a - g r )T ig a \. It can be noticed that in this weakly damped limit 



the integrands goes like 1/rj and the energy interval AVl — a < M < A give the main 
contribution to the current. 

Finally, when summing up all ac components to obtain the total dissipative and non- 
dissipative parts, the current densities become singular at the condition (p(u) = 0. This 
condition is satisfied for uj = us = iA^l — asin 2 (0/2) (i.e. at the bound states energy 
levels) leading to 



_ , e 2 a 2 A 4 r°o 2 « w 



(5uj sin 2 (fi(co) 
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(u - \u> s \ - iri){uj + \u s \ - irj) 
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Is(<f>) = — —a A 2 sin <p / np(o;)/m < 

h J-oo { (u — 
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\u s \ - irj){u + \u s \ - irj) J ' 

These integrals can be straightforwardly evaluated. Eq. (Cll) gives the expression for 
the phase-dependent linear conductance given in section IV (Eq. (34)), while Eq. (C12) 
yields the well known expression for the supercurrent in a single mode SQPC [13|,37[ 



_ / i \ eA 2 a sin(0) 

w = 9ft i ! i tanh 

2ft |^<?(0)| 
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Non-linear regime (r/ <C eV) 

We first rewrite Eqs. (A13-A15) for the current components as 
Jo = ^ /" dwtanh f ^ ) 2 ^ e {rr(a,7j(flg - 

ft V ^ / n=odd>0 

(ft ~ 9o)Tn9 a - n )} 



^ = ^r^tanh(^) £ i?e{Tr(a 
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where a rigid shift of nuo/2 in the energy arguments of the different T n with respect to the 
ones appearing in Eq. (A13-A15) has been introduced. 

In the limit eV — > 0, the solution of Eq. (A9) is A+ = X[u + (n + 2)a>o/2], where A (a;) 
satisfies the quadratic Eq. (C4). The coefficients T n (n > 0) can then be generated starting 
form T\ and using the transfer matrix z + (u). These quantities are obtained from Eqs. (A8) 
and (A10) making use of the eV —>■ solution for A+. 

The resulting expressions simplify considerably in the ballistic limit where A+ = i/t. For 
energies inside the gap one obtains 



f n {u) 



t 



n-l 



Z 3=1 



(C15) 
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where a n = arccos[(u; + neV)/A]. As discussed in Ref. ||, when written as in Eq. (C14), 
the main contribution to the current in this limit comes from a small energy range around 
the gap edges. Evaluation of these integrals leads directly to Eq. (35). 
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FIGURES 



FIG. 1. Schematic representations of the two different situations discussed in Section II. 

FIG. 2. The dc current-voltage characteristic of a S-S contact for different values of the normal 
transmission at zero temperature. 

FIG. 3. The first three ac components of the dissipative current for different values of the 
normal transmission at zero temperature. 

FIG. 4. Same as Fig. 4 for the non-dissipative current. 

FIG. 5. Current density for the dc component lo within the linear regime discussed in section 
IV. C. Cases (a), (b) and (c) correspond to transmission values a = 1, 0.65 and 0.04 respectively. 
In all cases the full line corresponds to 77/A = 1/10, the dotted line to 77/A = 1/25 and the broken 
line to 77/A = 1/100. The thermal factor sech 2 (/3u;/2) (see Eq. (C8)) has been extracted from the 
current density. 

FIG. 6. Behavior of the first non-dissipative (a) and dissipative (b) ac components in the very 
small voltage range close to ballistic conditions. These results have been obtained for negligible 77. 
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